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Abstract 

In mixture model-based clustering applications, it is common to fit several models from a family and 
report clustering results from the 'best' one. Selection of this best model is a difficult and consequential 
problem and criteria commonly used include the Bayesian information criterion, the Akaike information 
criterion, and the integrated completed likelihood. We propose an alternative to the selection of a best 
model, instead averaging the clustering results of several models. In the course of model averaging, 
the top few models often have different numbers of mixture components and so merging components 
is necessary. The effectiveness of our model-based clustering averaging approach is illustrated using a 
family of Gaussian mixture models on simulated and real data. This paper is perhaps the first step in 
a departure from the 'single best model' paradigm that currently dominates the model-based clustering 
literature. 



1 Introduction 



Model-based clustering is an idiom that describes the application of a mixture model, or any model, for 
clustering. Dating at least as far back as |Wolfe| (1963), interest in model-based clustering is increasing 
steadily in application areas like food authenticity, social networks, and microarray gene expression analysis 



(see 


Yeung et al. |2001| Wehrens et al. 


2004 


Krivitsky et al. 


2009 


McNicholas and Murphy 


2010 


Murphy 


et al. 


2010 


for examples). A random vector X arises from a parametric finite mixture distribution if its 



proportions (tt^ > 0, ^^^-^ tt^ — 1), and •& = (tti, . . . , ttg, 0i, . . . , 0g) denotes the parameters. Gaussian 
mixtures have dominated the model-based clustering literature until very recently. The likelihood for p- 
dimensional xi , . . . , x„ from a Gaussian mixture model is 



(1) 



where (/)(xi|/Xg, Sg) is the density of a multivariate Gaussian distribution with mean \ig and covariance Sg, 
and '& once again denotes the model parameters. 

The density in Equation [l] has a total of (G — 1) + Gp + Gp(p+ 1)/2 free parameters and so it is practically 
necessary to introduce parsimony. A total of Gp(p-)-l)/2 of these parameters are in the component covariances 
and so imposing the isotropic constraint on the component covariances, i.e.. Tig = Uglp, is a very simple 
way to reduce the number of parameters from quadratic to linear in p. Of course, this constraint will not 
be practical for most applications and so less restrictive constraints need to be considered. Several such 
approaches have been tried, usually based on imposing constraints upon a decomposed covariance structure. 
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The most famous such approach is based on an eigen-decomposition (Banfield and Raftery 

Sg = AgDgAgD^, WhCrB Dg 

elements proportional to its eigenvalues, and Xg 



1993), so that 



is the orthogonal matrix of eigenvectors of Sg, 
is a scalar. 



Ag is a diagonal matrix with 



Celeux and Govaert 



( 1995 1 impose constraints 

on Ag, Dg, and Ag to obtain a family of 14 mixture models. The famous MCLUST family of models 
(T able [T|) is the subset of these 14 that is supported by the mclust package (Fraley and Raftery 20061 for 



R (R Development Core Team 2012). Parameter estimation for each member of the family is carried out 



using the expectation-maximization (EM) algorithm (Dempster et al. 1977); extensive details on the EM 
algorithm and applications to mixture models are given by McLachlan and Krishnan ( 2008 ) . 



Table 1: Models available in the mclust package, along with number of covariance parameters in each case. 





Dist. 


Volume 


Shape 


Orient. 




Free cov. paras. 




EII 


Spherical 


Equal 


Equal 


NA 


AI 


1 




VII 


Spherical 


Variable 


Equal 


NA 


Agl 


G 




EEI 


Diagonal 


Equal 


Equal 


Cord. Axes 


AA 


P 




VEI 


Diagonal 


Variable 


Equal 


Cord. Axes 


AgA 


P + G-1 




EVI 


Diagonal 


Equal 


Variable 


Cord. Axes 


AAg 


Gp-G + 1 




VVI 


Diagonal 


Variable 


Variable 


Cord. Axes 


AgAg 


Gp 




EEE 


Ellipsoidal 


Equal 


Equal 


Equal 


ADAA' 


p{p+l)/2 




EEV 


Ellipsoidal 


Equal 


Equal 


Variable 


ADgAA; 


Gpip + l)/2-{G- 


-l)p 


VEV 


Ellipsoidal 


Variable 


Equal 


Variable 


AgDgAA; 


Gp(p + l)/2-(G- 




VVV 


Ellipsoidal 


Variable 


Variable 


Variable 




Gp(p + l)/2 





A typical application of the MCLUST family of models consists of running each of the ten models 
(Table [l]) for a range of values of G; the default in mclust is G = 1, . . . , 9, which results in 90 models being 
fitted. As is typical with families of mixture models, the best of these 90 models will be selected using 
some criterion and the associated classi fications report ed. The most popular criterion for this purpose is 
the Bayesian information criterion (BIC; 



Schwarz 



1978), BIC = 2?(x, i9) — plogn, where i? is the maximum 



likelihood estimate of i?, Z(x, t?) is the maximized log- likelihood, and p is the number of free parameters. 



Leroux (19921 and Keribin (2000) give theoretical results that, under certain regularity conditions, support 
the use of BIC for choosing the number of components in a mixture model. In addition, in a range 



of 



applications of model-based clustering, model selection based on the BIC has given good results (cf. Kass 



and Wasserman 


1995 


Dasgupta and Raftery 


1998 


Fraley and Raftery 


1998 



None the less, it is well known that the model selected based on BIC does not necessarily give the best 
predicted classifications or the most accurate number of components. 

Several alternatives to the BIC have been proposed. The Akaike information criterion (AIC), which 



was considered in the mixture modelling context by Bozdogan and Sclove ( 1984 ) , takes the form AIC = 
2Z(x, •&) — 2/9. However, it has been observed that the AIC is order inconsistent and tends to over-fit models 
(see Koehler and Murphree 1988 for example). 'Bicrnacki et al. (2000) proposed the integrated completed 
likelihood (ICL), which is approximated by penalizing the BIC using the estimated mean entropy. We have 



ICL « BIC + Eti MAP {z,g} log where 



1 if Xi belongs to the gth component and Zi, 







otherwise, and MAP denotes the maximum a posteriori probability such that MAP {z^g} = 1 if maxg {z^g} 
occurs in the gth component and MAP {zig} = otherwise. The ICL is an attempt to shift the focus 
from selection of the correct number of components G to selection of the right number of clusters. Another 
way to do this is to use the BIC, or another criterion, for model selection and then to merge the resulting 



Gaussian components to give clusters (cf. [Hennig 20101. |Li| ( j2005) proposed fitting a multilayer mixture 
(i.e., a mixture of Gaussian mixtures), assuming that the number of clusters k is known in advance, and 
then applies fc-means clustering to the G component means, where G is the estimated number of Gaussian 



mixture components based on the BIC. However, the method suggested by Li ( 2005, ) has the drawback that 
it is assumes a priori knowledge of the number of clusters and it does not meaningfully take account of 



2 



cluster shapes. Baudry et al. (2008) proposed fitting a Gaussian mixture model to data, then selecting the 



total number of Gaussian mixture components using BIG, and combining them hierarchically according to 
an entropy criterion. 

The goal of this paper is to depart from the paradigm wherein a single best model is selected. Instead, 
we report clustering results based on an averaging of the top few models. The remainder of this paper is 
outlined as follows. In Section [2] we describe our model averaging schemata, before illustrating our approach 
on real data (Section [3|. The paper concludes with discussion and suggestions for future work (Section |4]). 



2 Methodology 
2.1 Overview 

As mentioned earlier, it is common to fit many mixture models within a family and then report clustering 
results based only on the best one. Such criterion-based methods of model selection have the general feature 
that the larger the value of the criterion, the stronger the 'evidence' for the model (i.e., covariance structure 
and number of components, as well as number of latent variables where relevant). We argue that criterion- 
based approaches are not reasonable because they are effectively throwing away all but the best model. 
One may argue, inferentially, that it has to matter that 90 (say) models were fitted when reporting results 
based on just one. Practically, criterion-based approaches can be considered unreliable when the difference 
between the largest value of a criterion is close to other values. To overcome this problem, we consider a 
model averaging approach. More specifically, we consider an averaging of the top few mixture models with 
applications in clustering and classification. Herein, we consider cases where the number of models to be 
averaged is allowed to vary and cases where it is fixed a priori. 



2.2 Top Model Selection 

Herein, a top models set S is defined and will contain the 'top models' to be averaged. The selection 
procedure for the top models is the same for each of the three criteria, and in what follows here the BIG is 
utilized for illustration. Prior to the selection of the elements of 5, the top models are ordered according to 
their BIG values, e.g., denote the largest BIG value by BIG [1], and the corresponding model is the 1st model; 
the second largest BIG value is BIG [2], and the model is the 2nd model; likewise for the other models. Two 
scenarios are now considered: selection of the 'top M' and 'top 5' models, respectively. 



Top M Model Selection 

1. Step 1: Select the 1st model, and put it into S. Set fc = 2, 

2. Step 2: Select the fcth model, and calculate Diff= BIG [1] - BIG [k]. 

3. Step 3: If Diff< 10, put the fcth model into S, set fc ^ fc + 1, and return to Step 2. Else terminate. 
The resulting set S contains M models, i.e., the top M models. 



Top 5 Model Selection 

1. Select the top 5 models (corresponding to the first five highest value of the BIG), and put them into 
set S. 

The resulting set S contains 5 models, i.e., the top 5 models. 
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2.3 Merging Mixture Components 

If models within the set S have different numbers of components, we will need to merge the Gaussian mixture 
components. In this section, we introduce a procedure for merging a larger number of mixture components 
into a smaller number of clusters. For the model with larger number of components, e.g. G, and the new 
model after the merging process has smaller number of clusters, e.g. H, where H < G, the density takes 
the form /(x) = E.ti ^;/;(x) = Eg=i 

7rg(/)(x I /Xg, Sg), where tt* is the sum of TTg of the Gaussian mixture 
components assigned to f* , and f* (x) is from a single or a mixture of component densities from the original 
model, such that each original Gaussian component appears in exactly one of /j'(x), . . . , /^(x). Suppose, for 
example, that G = 3 and H ~ 2. Then /j* (x) could be a mixture of the first two Gaussian components and 
(x) would equal the third one. 



We use the adjusted Rand index (ARI; Rand 1971 Hubert and Arable 1985) in our merging criteria. 



The ARI is the Rand index corrected for chance. The Rand index compares two data partitions and is given 
by the number of pair agreements divided by the total number of pairs. The correction leading to the ARI 
is performed to account for the fact that if classification is performed randomly then some cases will be 
correctly classified by chance. 



The same notation from Section 2.2 is used here, where BIG [1] refers to the largest BIG value, and the 
corresponding model is referred as the 1st model, etc. We merge components to make the other models in 
S comparable to a 'reference model'. Two cases are considered. 



In Gase I, the 1st model has G components and is the reference model. All remaining models in set S 
are to be merged to G clusters. 

In Gase II, the ith {i ^ 1) model has the least number of components G and is the reference model. 
The remaining models (including the 1st model) are merged to have G components. 



Note that we are assuming that the BIG will only overestimate the number of components; therefore, we 
ignore models with less components than the reference model. 

The merging process is best illustrated with an example. Consider merging 7 components into 4 compo- 
nents, where 7 is the number of components from a model in S and 4 is the least number of components in S. 
Suppose we work under Case II, so the reference model has G = 4 components. The partition corresponding 
to the reference model is used as the underlying true classification in the merging process and is called the 
'reference partition'. In what follows here, to avoid confusion, we use {1, 2, . . . , 7} to denote components of 
the seven-component model and {a, 6, c, d} to denote components in the G = 4 component reference model. 
The following steps illustrate our merging strategy. 

1. A combination matrix A of size 35 x 4 is generated, note that (^) — 35, where each row represents 
a partial clustering, e.g., a sample row a22- = (2,3,4,6) means that group 2 goes into new group a, 
group 3 goes into new group 5, group 4 goes into new group c, and group 6 goes into new group d. 

2. For each row in A, we must deal with the remaining components; e.g., {1, 5, 7} in the case of row 22. A 
permutation matrix B (64 x 3) is created to include all possibilities. One sample row is hj. — (a, a, d), 
which puts component 1 into new component a, component 5 into new component a, and component 7 
into new component d, so the model after merging has the new groups {a, b, c, d} — {1U2U5, 3, 4, 6U7}. 

3. The ARIs between the reference partition and all partitions arising from the model after merging are 
calculated and recorded in the 35 x 64 matrix C. For example, in step 1, a22- = (2,3,4,6) is the 22nd 
row in matrix A, in step 2, hj. = {a,a,d) is the 4th row in matrix B, and the ARI value between 
reference partition and this model after merging is stored at the 22th row and the 4th column of the 
matrix C. 

Once we have every element of the matrix C, the best merging combination is chosen to correspond to the 
largest ARI value. 
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2.4 Top Model Averaging 



After merging components where needed so that each model within S has the same number, and discarding 
those with less components than the reference model (if applicable), it remains to find the estimated con- 
ditional probability for the model after merging and then to average these conditional probabilities across 
the models in S. Denote by Zf^ the conditional posterior probability that arises from component g with 
respect to the G-component solution. If (only) components g' and g" from the G-component solution are 
combined or merged into one cluster, then the Zig remain the same except for g' and g", and the new 
component g' U g" has the conditional probability Zig>ug" = Zf^, + Z^„. 

Sticking with the notation from Section |2.2[ model averaging is a process to compute the weighted 
summation of the models in the set S, with Zavo = ^W-^Mi where Zi^^e is the averaged model conditional 
probability matrix, w [i] is the weight of the ith model defined later, and Z [i] is the conditional probability 
matrix for the ith model. The purpose of the w [i] is to weight the importance of each model. Essentially, 
the higher the value of the criterion (BIG, ICL, or AIC), the stronger the evidence for the model and the 
number of components, and the larger this weight should be. In this regard, w [i] is defined as 



where v [i] is the value of the criterion for the ith model. Note that higher v [i] values produce larger w 



3 Data Analyses 
3.1 Italian Wine 



Forina et al. ( 1986 1 present data on chemical and physical properties of 178 samples of three varieties (Barolo, 



Barbera, and Grignolino) of wine from the Piedmont region of Italy. A data set containing 27 physical and 



chemical properties is available in the pgmin package ( McNicholas et al. 2011 1 for R. We ran mclust on these 



data, under the default settings. The VVI model with G — 3 components was the best model in terms of 
BIG (-12103.76) and the EVI model had almost identical BIG (-12103.82). The IGL chose the same top 
two models, albeit in the opposite order, and the AIG selected a G = 8 component model. We then applied 
our top 5 and top M averaging procedures (merging components where needed). The results for the top 5 
averaging approach (Table [2]) show that averaging led to much better clustering performance (taking the 
wine regions as the correct classes). This is true for the BIG (ARI of 0.914 versus 0.895), the IGL (ARI of 
0.914 versus 0.830), and the AIG (0.722 versus 0.386). The results for the top M approach (Table [s]) are 
not as impressive, with averaging performing as well as or worse then selecting the best model. Note that 
M was small for all three criteria: M = 2 for BIG and IGL, and Af = 1 for AIG. 

3.2 Italian Olive Oil 



Forina and Tiscornia ( 1982 ) present data on the percentage composition of eight fatty acids found by lipid 
fraction of 572 olive oils from three regions of Italy (Southern Italy, Sardinia, and Northern Italy). These 
three regions can be further broken down into nine areas. Again, these data are available in the pgmm package. 
As before, we ran mclust on these data using the default settings. The BIG and the IGL both selected a 
VVV model with G = 5 components. The AIG selected a VVV model with G = 9 components. We then 
applied our top 5 and top M averaging procedures, merging components where needed. The results for the 
top 5 averaging approach (Table [4| show that averaging led to superior clustering performance no matter 
whether the regions or the areas were taken as the true classes. The top M approach could not be applied 
to these data because the difference between the best and subsequent models was more than 10 under each 
criterion. 
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Tabic 2: Results for top 5 model selection, merging, and averaging for the wine data. 





Model Selection 


Merging 


Average 

VV Cl^lit 


ARI 


BIC 


Ist(VVI) -12103.76 3 
2nd(EVI) -12103.82 3 
3rd(VEI) -12145.07 7 
4th(EEI) -12191.58 5 
5tVil'F,ET1 -191 93 69 6 


no need 
no need 

{1 U 2, 3 U 4 U 5, 6 U 7} no need 

{1,2U3,4U5} 
1 1 1 1 9 3 114 1 1 6]- 


0.2007217 
0.2007208 
0.200039 
0.199276 
1QQ9495 


Best: 
0.8950872 

Averaged: 

Q1 35014 


ICL 


Ist(EVI) -12104.4 3 
2nd(VVI) -12104.45 3 
3rd(VEI) -12149.27 7 
4th(EEI) -12193.59 5 
5th(EEI) -12198.58 6 


no need 
no need 

{1 U 2, 3 U 4 U 5, 6 U 7} no need 

{1,2U3,4U5} 
{1U2,3U4,5U6} 


0.2007521 
0.2007513 
0.2000107 

0.1992837 
0.1992022 


Best: 
0.8301066 

Averaged: 
0.9135014 


AIC 


Ist(EVI) -8862.6 8 

2nd(EVI) -9152.212 7 
3rd(VEI) -9643.861 7 
4th(EEI) -10052.09 6 
5th(VEV) -10432.34 5 


{1U2U3,4,5U6,7, 8} 
{1U2,3,4U5,6,7} 
no need {1 U 2, 3,4 U 5, 6, 7} 
{1U2,3,4,5,6} 
no need 


0.2165185 

0.209667 
0.198978 
0.1908972 
0.1839393 


Best: 
0.3860865 

Averaged: 
0.7215702 



Table 3: Results from top M model selection, merging, and averaging for the wine data. 





Model Selection 

Models Criteria G 


Merging 
Cases I & II 


Average 
Weight 


ARI 


BIC 


Ist(VVI) -12103.8 3 
2nd(EVI) -12103.8 3 


no need 


0.5000012 
0.4999988 


Best: 0.8950872 
Averaged: 0.8461782 


ICL 


Ist(EVI) -12104.4 3 
2nd(VVI) -12104.5 3 


no need 


0.500001 
0.499999 


Best: 0.8301066 
Averaged: 0.8301066 


AIC 


Ist(EEV) -8862.6 8 


no need 


no need 


Best: 0.3860865 



Table 4: Results for top 5 model selection, merging, and averaging for the olive data. 





Model Selection 
Models Criteria G 


Merging 

Case I Case II 


Average 
Weight 


ARI 

Regions Areas 


BIC 


Ist(VVV) -5096.72 5 
2nd(WV) -5166.82 9 

3rd(VVV) -5225.58 8 
4th(VVV) -5228.28 6 
5th(VVV) -5285.89 7 


no need 
{1U2,3U4,5U6,7,8U9} 
{1 U 8, 2, 3 U 4, 5 U 6, 7} no need 
{1,2U3,4,5,6} 
{1,2U3,4U5,6,7} 


0.2040479 
0.2012789 

0.1990155 
0.1989128 
0.1967449 


Best: Best: 
0.5969466 0.776304 

Averaged: Averaged: 
0.6057673 0.8056072 


ICL 


Ist(VVV) -5098.71 5 
2nd(WV) -5177.5 9 

3rd(VVV) -5235.94 8 
4th(VVV) -5250.09 6 
5th(VVV) -5307.68 7 


no need 
{1U2,3U4,5U6,7,8U9} 
{1U8,2,3U4,5U6,7} no need 
{1,2U3,4,5,6} 
{1,2U3,4U5,6,7} 


0.2044837 
0.2013719 

0.1991243 
0.1985875 
0.1964326 


Best: Best: 
0.5969466 0.776304 

Averaged: Averaged: 
0.6057673 0.8056072 


AIC 


Ist(VVV) -3409.77 9 
2nd(WV) -3664.24 8 
3rd(VVV) -3920.26 7 
4th(VVV) -4058.36 6 
5th(VVV) -4122.5 5 


{1U2,3U4,5U6,7,8U9} 
{1U8,2,3U4,5U6,7} 
no need {1, 2 U 3, 4 U 5, 6, 7} 
{1,2U3,4,5,6} 
no need 


0.2238324 
0.2082877 
0.1946852 
0.1880603 
0.1851345 


Best: Best: 
0.3437476 0.6600636 

Averaged: Averaged: 
0.6057673 0.8056072 
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3.3 Simulated Data 



Data were simulated from a mixture of G = 4 multivariate Gaussian distributions using the genRandomClust () 
function from the R package clusterGeneration (Qiu and Joe. 2012 1. The genRcindoniClust ( ) function 



generates random clusters based on the method proposed by Qiu and Joe ( 2006[ ) . We generate p — 3 vari- 
ables for n = 562 observations and we choose to generate substantially overlapping clusters (Figure [T]). This 
is a very difficult clustering problem and perfect classification is not expected. We ran mclust on these data, 
under the default settings. The EII model with G = 5 was selected by the BIG, the ICL selected a G = 4 
component EVI model, and the AIG chose a G = 9 component EVI model. 




Figure 1: Pairs plot of the G = 4 component, p = 3 variable data generated by genRandomClust () . 

We then applied our top 5 and top M averaging procedures, merging components where needed. The 
results for the top 5 averaging approach (Table [s]) show that averaging led to far superior clustering perfor- 
mance when the BIG is used, taking the ARI from 0.714 for the single best model to 0.823 (Gase I) or 0.878 
(Gase II). Top 5 merging had no effect for the IGL but did increase the ARI when the AIG was used (from 
0.535 to 0.593). The results for the top M approach (Table |3]) are not applicable for the BIG or the IGL but 
had a dramatic effect when used with the AIG. A total of 17 models had an AIG value within 10 of the best 
model and averaging these 18 clustering models, with Gase II merging, resulted in a remarkably improved 
ARI (0.881 versus 0.535). 

4 Summary and Discussion 

This paper is the first step in what will perhaps ultimately be a departure from the 'single best model' 
paradigm that has heretofore dominated the model-based clustering literature. The averaging approaches 
used are not technically complicated but we show that they can produce extraordinary improvement in 
clustering performance over choosing a best model. At times, our merging approach requires that components 
are merged and we introduced methods to do this. We explored two options for determining the number of 
components to merge over: we can choose the number of components attached to the best model (Case I) 
or we can chose the model in S with the lowest number of components (Case II). 

The approach introduced herein is just the beginning and the are tremendous opportunities to build upon 
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Table 5: Results for top 5 model selection, merging, and averaging for the simulated data. 





Model Selection 




Mer 


ging 


Average 


ARI 




Models 


Criteria 


G 


Case I 


Case II 


Weight 






Ist(EII) 


-9753.587 


5 


no need 


{1,2U3,4,5} 


0.2002342 


Best: 0.7135163 




2nd(EVI) 


-9763.885 


4 


no need 


no need 


0.200023 


Avrg. (Case I): 


BIG 


3rd(EEE) 


-9767.062 


6 


{1,2,3, 4,5U6} 


{1,2U4,3,5U 6} 


0.1999579 


0.8227513 




4th(EII) 


-9768.561 


7 


{1,2U3,4,5,6U7} 


{1,2U3U5,4,6U7} 


0.1999272 


Avrg. (Case II): 




5th(EEI) 


-9771.959 


7 


{1,2U3,4,5,6U7} 


{1,2U3U5,4,6U7} 


0.1998577 


0.8777635 




Ist(EVI) 


-9795.119 


4 


no need 




0.2005361 


Best: 




2nd(EEV) 


-9812.766 


4 


no need 




0.2001755 


0.8760625 


ICL 


3rd(EEE) 


-9822.785 


6 


{1,2U4,3,5U6} 


no need 


0.1999713 






4th (E VI) 


-9835.585 


5 


{1,2,3,4U5} 




0.1997111 


Avrg. (Case I): 




5th(EII) 


-9840.764 


5 


{1,2,3,4U5} 




0.199606 


0.8760625 




Ist(EVI) 


-9638.797 


9 




{1, 2, 3, 4, 5, 6 U 7, 8, 9} 


0.2000247 


Best: 




2nd(EII) 


-9639.507 


8 




no need 


0.20001 


0.5347199 


AIC 


3rd(EVI) 


-9640.244 


8 


no need 


no need 


0.1999947 






4th(EEI) 


-9640.36 


8 




no need 


0.1999923 


Avrg. (Case I): 




5th(VEI) 


-9641.026 


8 




no need 


0.1999784 


0.5931062 



Table 6: Top M model selection, merge and average for the simulated data. 





Model Selection 






Merge 


Average 


ARI 




Models 


Criteria 


Clusters 


Case 1 


Case 2 


Weight 




BIC 


Ist(EII) 


-9753.587 


5 


no need 


no need 




Best:0.7135163 


ICL 


Ist(EVI) 


-9795.119 


4 


no need 


no need 




Best:0.8760625 




Ist(EVI) 


-9638.797 


9 




{1U9,2U3U5,4U8,6,7} 


0.05558314 






2nd(EII) 


-9639.507 


8 




{1,2U3U5,4U8,6,7} 


0.05557905 






3rd(EVI) 


-9640.244 


8 




{1,2U3U5,4U8,6,7} 


0.0555748 






4th(EEI) 


-9640.36 


8 




{1,2U3U5,4U8,6,7} 


0.05557413 


Best: 




5th(VEI) 


-9641.026 


8 




{1,2U3U5,4U8,6,7} 


0.05557029 


0.5347199 




6th(EEE) 


-9641.448 


6 




{1,2U4,3,5,6} 


0.05556786 






7th(EEI) 


-9642.014 


7 




{1,2U3U 5,4,6,7} 


0.0555646 






8th(EVI) 


-9642.498 


7 




{1,2U3U5,4,6,7} 


0.05556181 




AIC 


9th(EII) 


-9643.443 


9 


no need 


{1U9,2U3U5,4U8,6,7} 


0.05555636 






lOth(EEE) 


-9643.724 


7 




{1,2U3U 5,4,6,7} 


0.05555474 


Avrg. (Case II): 




llth(VEI) 


-9644.026 


7 




{1,2U3U5,4,6,7} 


0.055553 


0.8811323 




12th(EEE) 


-9644.5 


8 




{1,2U3U5,4U8,6,7} 


0.05555028 






13th(VII) 


-9645.838 


9 




{1U9,2U3U5,4U8,6,7} 


0.05554257 






14th(EEI) 


-9646.4 


9 




{1U9,2U3U5,4U8,6,7} 


0.05553934 






15th(EEE) 


-9646.47 


9 




{1U9,2U3U5,4U8,6,7} 


0.05553893 






16th(EII) 


-9647.279 


7 




{1,2 U3U 5,4,6,7} 


0.05553427 






17th(VVV) 


-9648.39 


6 




{1,2U4,3,5,6} 


0.05552788 






18th(EVI) 


-9648.552 


5 




no need 


0.05552695 





it. We illustrated merging for three different criteria (BIC, ICL, and AIC) but the same approach could 
be used with other criteria (e.g., the LASSO-penalized BIC of Bhattacharya and McNicholas 2012 1. We 
used Gaussian mixture models but the same merging procedures could be applied to non-Gaussian mixtures. 
We focused on clustering applications, but our merging methods could also be applied for model-based 
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classification (cf. Dean et al. 2006 McNicholas 20101. The idea of a top 5 approach, while we have shown 



it can be effective, is a little ad hoc and so is the use of a criteria value difference of 10 in the top M 
approach. One could, for example populate S with all models with a criterion value within some proportion 
of the best model. The choice of weight is also flexible. The weight definition used herein achieved improved 
clustering performance but it needs to be tested extensively to better establish its effectiveness. This, along 
with everything else mentioned in this paragraph, will form the basis of future work. 

Some may criticize this approach because the final result is not an interpretable model in the sense that 
a Bayesian model averaging might produce an interpretable model. They may go further and argue that 
the main advantage of mixture model-based clustering is lost when one averages, as we do, placing value on 
classification accuracy ahead of model interpretability. Of course, as mentioned herein, the idea of reporting 
clustering results for one model as if the other 89, say, had not been run could also be considered questionable. 
Whatever about the philosophical debate, if the ultimate goal of clustering is to obtain the best possible 
estimated group memberships then we maintain that our results show that merging model-based clustering 
results is worth a further look. 
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